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Abstract 



, We present a detailed study of different exchange-correlation (xc) functionals 

^ \ in describing the dynamical properties of finite systems. For that purpose, 

\ we calculated the static polarizabilities, ionization potentials and optical ab- 

sorption spectrum of four small clusters, Na2, Na4, SiH4 and Si2H6, using 
CN . a real-space, real-time technique. The computed static polarizabilities and 

, ionization potentials seem to be in rather good agreement with the available 

\ experimental data, once the proper asymptotics of the potential are taken 

I into account. The same conclusion holds for the absorption spectra, although 

■ the xc kernels in use do not provide a sufficiently strong attractive interac- 
tion between electrons and holes, leading to spectra slightly shifted towards 
higher energies. This deficiency is traced back to the insufficient description 

\ of dynamical effects in the correlation functional. Furthermore, it is shown 

that the xc potential used to obtain the ground state is the key factor to get 

■ reasonable spectra, whereas the choice of the xc kernel just amounts to small, 
^ ■ although important, quantitative changes. 
O 

I. INTRODUCTION 

■ Density Functional Theory (DFT) became, in the two last decades, the method of 
election for the ah initio calculation of material properties. Systems with a couple thousand 
atoms are now routinely investigated, and the calculated energies, geometries, etc., very often 
agree spectacularly with the experimental data. This success, however, was not immediate: 
More than twenty years of research were necessary to obtain an exchange-correlation (xc) 
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energy functional which was precise enough to satisfy the quantum-chemistry needs. The 
accurate and sophisticated Generahzed Gradient Approximations (GGA) that we now have 
at our disposal are indeed the result of a long history of attempts, tests and failures. 

Despite these remarkable achievements, there are some quantities which are beyond the 
reach of the conventional DFT theory. In this Article we will be concerned by one of such 
properties, namely electronic excitation energies and optical spectra. Several extensions 
of the original framework were put forward to obtain excited-state properties. From all 
different approaches, Time-Dependent DFT (TDDFT) p[ emerged the last years as the 
main DFT formalism for these calculations. One of the main advantages of the TDDFT 
formulation is that it allows one to deal with problems beyond the perturbative regime as, 
for example, the response of atoms to ultra-short and intense laser beams. In this way, 
by solving the time-dependent Kohn-Sham (KS) equations, we can get, e.g., information 
concerning the harmonic spectrum or the ionization yields. This type of approach is very 
important in condensed matter science where femtosecond laser pulses are used to monitor 
the dynamics of electrons in a solid. Furthermore, with larger computational power, we may 
deal with non-adiabatic couplings between the electronic and the ionic degrees of freedom 
in presence of these high electromagnetic fields (we may observe phonon assisted structural 
transformation induced by these external fields). 

Excitation energies can be obtained from TDDFT either from the position of the poles 
of the KS linear response function or from the time-dependent density. In the second 
case, we will be using an approximated xc time-dependent potential, fxc(r,t), while in the 
former, the key quantity is the xc kernel, /xc(r, r', t, t'), defined by the functional derivative: 



As expected, the xc kernel, due to its complicated structure, is much harder to model than 
the local time-dependent xc potential. Furthermore, both calculations usually start from 
the KS ground-state of the system, which depends on the approximation used for the static 
Wxc(r)- 

One of the most widely used approximation to time- dependent phenomena is the adi- 
abatic local density approximation (ALDA). Although this functional is constructed using 
the xc correlation energy of the (constant density) electron gas, it yields, like the LDA 
does in the static case, rather accurate results for systems with rapidly varying densities, as 
atoms, surfaces or clusters For example, the photo-absorption spectrum of rare-atom 

gases was computed in Ref. and the agreement with experiment is remarkably good. 
Results of similar quality have been achieved for the photo-response of small metallic and 
semiconducting clusters |0. For the metallic clusters it was shown that the inclusion of xc 
effects in the dielectric response was important to get the correct red-shift, as compared to 
experiments, while for the semiconducting (silicon) clusters it was found that the spectra 
of different isomers was sufficiently different to distinguish between them. Surface and con- 
finement effects were responsible for the appearance of absorption in the optical range for 
the two silicon clusters, Si4 and Sig. Similar results were found for C20 clusters 0. 

The purpose of this work is twofold: To address the impact of a good xc-potential in 
the optical spectrum of small clusters and to estimate how relevant the /xc kernel is in this 
respect. We note that in Ref. the quality of different functionals was tested for light 
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atoms, and in Ref. |T0| for the calculation of the correlation energy of the homogeneous 
electron gas. In the first of these works, which included calculations using the exact static 
xc potential, the correct description of the xc-potential used to obtain the ground state 
seemed to be the dominant factor, whereas the kernel played a marginal role. However this 
was not the case when looking at the total correlation energy of the homogeneous electron 
gas. In particular, the nonzero spatial range of fxc{r,r',uj) could not be neglected, whereas 
the frequency dependence appeared to be less important. Our work clearly demonstrates 
that the optical spectrum of small systems is not only determined by the static xc potential, 
but also the kernel plays a role. 

The rest of the Article is structured in the following way: In Section we give a 
brief overview of the xc functional zoo used in static KS calculations; We then proceed by 
explaining how these functionals can be changed in order to be used in TDDFT calculations; 
In Section |I^ we give some details on the techniques used in our calculations, and present 
our results for some selected sodium and silicon clusters; We finally conclude and give some 
remarks on the quality of the functionals tested. 



II. XC FUNCTIONALS FOR STATIC DFT 

The first (and simplest) approximation for the xc functional, the Local Density Ap- 
proximation (LDA) introduced by Kohn and Sham in 1965 Jill, yielded remarkably good 
results for such a modest effort, but the method didn't prove accurate enough to be used in 
theoretical chemistry. The LDA xc energy functional can be written as 

E^^'^ = Jd\n{r)e^^{n{r)), (2) 

where n(r) is the electronic density at point r, and ^^c'^in) stands for the xc energy density 
of the homogeneous electron gas with density n. This function is calculated very accurately 
using Monte-Carlo techniques |jl2|, and then fitted to some simple analytical function [jl3|. 



The functional (0) has several shortcomings, among which we point out the following: 

• It neglects non-local effects, i.e. the LDA xc energy density at point r only depends 
on the density at that point. We therefore should not expect this functional to work 
in cases where the density has very strong spatial variations. 

• The exchange part of the functional does not cancel exactly the self-energy part of 
the Hartree term. This leads to a wrong asymptotic behavior of the xc potential 
for finite systems (it goes exponentially to zero, instead as — e^/r). Properties that 
strongly depend on this asymptotic behavior, like the ionization potential of atoms 
and molecules, come out with very large errors. Also, there are no Rydberg states 
within the LDA, and negative ions usually do not bind, which renders the calculation 
of electron-affinities impossible. 

• LDA usually overbinds, giving too short bond-lengths, etc. 

The next generation of functionals included the so-called Generalized Gradient Approx- 
imated (GGA) functionals []T1|-|T^. They can be written as: 
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GGA 



n{v) e^^^{n{r), Vn(r)), (3) 



where Vn(r) is the gradient of the density at the point r. e^^^ is usually some analytic 
function with some free parameters that are either fitted to experiment, or determined 
by some exact sum-rules. GGAs solve some of the problems present in the LDA, and 
in some cases yield results with high enough precision to be used as a tool in quantum 
chemistry calculations. Recently, a new class of functionals generalizing the GGAs has 
been proposed. These so-called Meta-GGAs [|l7l, depend explicitly not only on the density 
and its gradient, but also on the kinetic energy density^. This extra dependence adds 
more flexibility, and allows better approximations to the exact xc functional to be built. 
Unfortunately, all the GGAs and MGGAs proposed till now suffer, to different degrees, 
from the same self-interaction problem as the original LDA. A way to circumvent it was 
proposed by van Leeuwen and Baerends in 1994 flSl: They applied Becke's construction 



not to the derivation of the exchange energy functional (like in Becke's original work), but 
to the modeling of the xc potential directly. By imposing the correct asymptotic behavior 
to the potential, they were able to get much better ionization potentials (and eigenvalues in 
general) . 

To complete the zoo of energy functionals, we still have to refer to orbital functionals, in 
the so-called Optimized Effective Potential (OEP), or Optimized Potential Method (0PM) 
T9| . This third generation of functionals is written explicitly in terms of the Kohn-Sham 



orbitals (being nevertheless implicit functionals of the density): 

ET = E^r[Mr)---Mr)]. (4) 
The xc potential is then calculated using twice the chain-rule for functional derivatives: 



J2 Jd'r' dh" 



5n{j) 



+ c.c. 



The first term of the right, f^^^^pj, can easily be obtained from the expression for 

E^^^ [({)i{v) ■ ■ ■ (^^{r)], while the second, ^^^^Ji^ can be calculated using first-order per- 
turbation theory. Finally, the remaining functional derivative can be identified with the 
inverse of the response function for non-interacting electrons, Xks('''''')- Rearranging the 
terms in Eq. (^, we obtain an integral equation for v^^^ . The solution of this integral 
equation is numerically very involved, and has only been achieved for systems with very 
high symmetry. An alternative is to perform an approximation first proposed by Krieger, 
Lee and lafrate (KLI), which transforms the hard task of solving an integral equation 
in three dimensional space into the simple one of solving a small set of linear equations. 



^ The MGGAs, due to the explicit dependence on the kinetic energy density are, in fact, not 
simple GGAs, but orbital functionals. 
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Two examples of these orbital functionals are the exact exchange (EXX), and self- 
interaction corrected (SIC) LDA functionals 0]. In the first, one employs the exact expres- 
sion for the exchange energy: 



E. 



EXX 



^EEJdrdr 



(6) 



o- i,j=l 



The SIC-LDA functional, originally proposed by Perdew and Zunger |13|, can be written as: 



N 



= E^f^ K(r), n,{r)] ^x°^ > 



c i=l 



N 



CF 1 = 1 



r — r 



(7) 



It is clear that the SIC-LDA obeys two of the features of the exact functional: It exactly 
cancels the self-interaction part of the Hartree energy and it vanishes for one electron sys- 
tems. The SIC-LDA functional is nevertheless ill-defined, for it is not invariant upon an 
unitary transformation of the Kohn-Sham wave-functions. 



III. XC FUNCTIONALS FOR TDDFT 



The simplest, and most commonly used, approximation to the xc functional in TDDFT 
is the adiabatic LDA (ALDA), in which the static LDA functional is used for the dynamical 
properties, but evaluated at the time- dependent density: 



t) 



(8) 



In the ALDA, the /xc kernel is a contact function in time and space: 



dv- 



hom/ 



n] 



dn 



(9) 



n=n(r,t) 



Following the same reasoning, it is straightforward to write adiabatic GGA potentials and 
/xc kernels. Unfortunately, the thresholds for the onset of absorption calculated either with 
adiabatic LDA or GGA functionals are typically below the observed ones (by several eV in 
the case of atoms). This is an intrinsic drawback of the ALDA/GGA because, in principle, 
TDDFT should yield the correct thresholds. This problem is once more related to the wrong 
asymptotic behavior of the effective Kohn-Sham potential that goes exponentially to zero 
instead of — e^/r for neutral systems. This is due, as already mentioned in the previous Sec- 
tion, to the insufficient correction of the self-interaction part of the Hartree potential. Also 
the xc kernel, /xc, suffers from a self-interaction error. A simple way to correct the asymp- 
totic part of the adiabatic potential is using the Adiabatic LB94 approximation. However, 
the high lying excitation energies calculated with this functional are usually overestimated 
for small molecules, and the LB94 does not do quite as well as the ALDA for the low-lying 



states 21 . 
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A generalization of the EXX potential to the time-domain is also possible [^. It starts 
from the perturbative expansion of the action functional and then uses the chain-rule for 
functional derivatives to derive an integral equation for the xc potential. A KLI-like approx- 
imation can then be used to simplify the task of solving this integral equation. Using the 
same procedure, one can derive an integral equation for /xc- However, for practical calcula- 
tions, it is desirable to devise a simple analytical expression for the xc kernel. To this end, 
it is common practice to further simplify the KLI potential by neglecting one of its terms 



/x 



. In the static exchange-only case, this leads to the so-called Slater approximation ^3 
then reads: 



EXX approx. 



rt, r't') = -6{t - t')5„ 



|r — r'l rio-(r)no-(r') 



(10) 



One should note that, due to the extra functional derivative in the definition of the xc kernel, 
it is much more complicated, in this case, to evaluate directly /xc than to time-propagate 
the KS states with t>xc. 

Finally, an adiabatic SIC-LDA can be easily written in the spirit of Eq. (|^. 

Other functionals designed to retain features which are important in time-dependent 
calculations have appeared in Quantum-Chemistry literature over the past years. A brief 
account of these can be found in Ref. 1211 • 



IV. METHOD AND RESULTS 

We used a real-space, real-time approach to solve the TD Kohn-Sham equations for 
Na2, Na4, silane (SiH4) and disilane (Si2H6). This time-evolution method is adequate to be 
combined with real-space calculations of the ground-state by means of the finite-difference 
pseudopotential method p4|,p5|, or with adaptative coordinates [Q. In both real- 



space discretization of the kinetic energy operator leads to sparse Hamiltonian matrices 
which do not need to be stored in memory, and are easy to handle. Furthermore, in order 
to propagate the states in time, we do not need to compute the complicated /xc kernel, but 
only the much simpler fxc potential. In this way, we are able to do time-dependent EXX 
calculations without having to perform Slater's approximation ||22| , p3| . 



The space was uniformly discretized, with mesh-spacing ranging from 0.4A for the sodium 
clusters to 0.25A for the silicon ones, and the points were contained inside a sphere with 
radius lOA for sodium and 7A for silicon. The ionic potential was modeled by a soft-core 
Troullier-Martins pseudopotential ^\ . By using these settings, a convergence of better than 



0.1 eV was reached for the eigenvalues. The time-evolution was performed with a modified 
Krank-Nicholson rule, with a time-step of 0.0025 h/eV, for a total evolution time, T^ax, of 
50 h/eV for the sodium clusters and 25 h/eV for the silicon clusters. With this parameters 
we can resolve excitation energies within AE = 27ih/Tmax that corresponds to tenths of 
eV. We added absorbing boundary conditions to improve the quality of the spectrum above 
the ionization threshold (avoiding artificial formation of standing waves in the spherical box 
used to confine the cluster). 

To calculate the dipole strength function (which is simply proportional to the absorption 
cross-section), we prepare our system in the ground-state, and then excite it with a delta 
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electric field, EoS(t). The dipole strength can then be related to the imaginary part of the 
dynamical polarizability: 

S{u) = ^uQa{u), (11) 

where h is Planck's constant, mg the electron's mass, and the dynamical polarizability, 

2 f 

a{uj) = —— dr z 6n{r,uj). (12) 

Eq J 

In the last expression, 6n{r, u) stands for the Fourier transform of n{r, t) — n(r, t = 0). The 
calculation of the dipole strength involves two important approximations: The choice of the 
static xc potential used to obtain the ground-state, and the time-dependent xc potential 
(directly related to the /xc kernel) used to propagate the state. The relative importance of 
both is still fairly obscure. Petersilka et al. linear-response results 0] for the Helium and 
Beryllium atoms indicates that the choice of the /xc kernel is less dramatic than the choice of 
a good static f xc potential. However, this is no longer true for the lowest excitation energies 
of Beryllium and if we look at the singlet-triplet splittings, where the effects arising from f xc 
cancel. The fairly accurate results of Lithium and Beryllium within the LDA were related 
to the quality of the ground state calculation for few electron systems as light atoms and 



the H2 molecule 28 



By performing the time evolution, in the limit t — 00, we have all the information about 
the linear response of the system, and should therefore be able to study the xc kernel without 
having to perform explicitly the functional derivative of the xc potential. Thus, we can study 
the effect of different /xc kernels for a fixed Vxc potential in the time-evolution method. To 
address this question, we propose the following method: First, we obtain the ground state 
of the system with the static functional v^J , and then perform the time evolution with the 
xc potential: 



vUr,t) = vg^{r,t = 0) + [4c^(r,t) - v^JJ{r,t = 0)J , (13) 

where v^^ is some other xc functional. In this way, we probe the importance of the /^^^ 
kernel, without having to spend more computational time in calculating this rather complex 
quantity. This method will be used for silane, for a varied selection of v^^ and v^J to 
quantify the dependence of the dynamical spectrum with the xc-potential and kernel. 

In the following, PZ stands for LDA |Tl[|, with ej?°™ taken from Monte-Carlo calculations 
and then parameterized by Perdew and Zunger PBE is the Perdew, Becke and 



Ernzerhof GGA functional LB94 represents the van Leeuwen and Baerends potential 
TE| ; SIC is the SIC-LDA functional |T3[; and, finally, EXX stands for exact-exchange . 



SIC and EXX are both treated within the KLI approximation |20]. "exp" will be used to 



denote experimental values. All energies are measured in eV unless otherwise stated. 



A. Sodium 

In Table | we show the ionization potentials (IP), obtained from the highest occupied KS 
eigenvalue, and HOMO-LUMO gap for the two sodium clusters. The IPs obtained either 
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by using the LDA or the GGA functionals are much smaller than the experimental value. 
This is a well known problem related to the wrong asymptotic behavior of the xc potentials. 
The situation does improve enormously if we correct for this deficiency, as can be seen from 
the LB94, SIC and EXX values. Although the calculated IPs differ by more than 2 eV, the 
HOMO-LUMO gap stays essentially constant, i.e, correcting for the asymptotic part of the 
exchange potential amounts to a constant shift in the eigenvalues, and not to an opening of 
the gap. 

All the calculated optical spectra for Na2 are quite similar (see Fig. |l|), regardless of 
the xc potential used, and exhibit three clear peaks in the 2-5 eV range (the third peak 
in the LB94 curve is almost completely smeared out). They all compare quite well with 
experiment, although the DFT peaks are all shifted towards higher energies. This shift can 
be understood in terms of the competition between the Coulomb repulsion in the electronic 
kernel and the electron/hole attraction from the xc part of the response (that is, Ae = 
€-c~^v+ < fv'Pcl |ri-r2| I'^^'V^e ^ +^xc)- In particular for a givcu valence-couductiou transition 
f — > c, the exchange correlation kernel within the simple ALDA introduced only a local and 
static attractive electron-hole interaction: 

= Jd\^:{T)Mr)^^^^:{r)^c{r). (14) 

This expression is more complicated for the other kernels, but it is clear that they do 
introduce an effective attractive interaction that it is not complete. We might have to 
recall that we have neglected dynamical effects and they may need to be included in the xc 
kernel to recover this minor effect. Also, we should keep in mind that temperature effects 
(vibrational motion) of the molecule may introduce a broadening of the spectrum as well as 
a shift of the peaks to lower frequencies ||7|,p9|. The functional which yields the best results 



for the dimer is, by a small margin, the EXX, while the strongest depart from experiment 
is the LB94 curve. We note that the 3rd peak in the LDA and PEE spectra lies in the 
continuum of states (above the ionization threshold). Thus it is more a resonance than a 
well defined bound transition as observed in experiments This deficiency is observed 



in all the clusters studied in the present work and it is, once more, mainly due to the wrong 
description of the asymptotic part of the exchange potential. 

For Na4 (Fig. |^) all DFT calculations yield very similar spectra (and similar to the 
many-body results based on a GW quasiparticle calculation and including the electron-hole 
interaction through the solution of the Bethe-Salpeter equation |^). The spectra consist of 



three well separated peaks in the 1.5-3.5 eV range, and a broader feature at around 4.5 eV. 
The comparison with the experimental values is quite good, although the peaks appear, once 
more, at higher energies. The best results were obtained using the LDA, the GGA and the 
EXX, while the LB94 spectrum showed the strongest deviation from experiment. Again, 
the error is larger for the high energy peaks, where transitions to conduction states close to 
the ionization threshold are involved. Due to the small overlap between these Rydberg-like 
states and the low-lying states, exchange processes are not relevant and the polarization con- 
tribution to correlation, although being generally weak becomes the dominant contribution 
to the renormalization of these single particle excitations. 
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Finally, in Table ||, we present the static polarizabilities of Na2 and Na4, obtained 
either through a finite field (static) calculation, or from the Fourier transform of the time- 
dependent density. For Na2 all functionals perform equally well, being the results smaller 
than the experimental value by around 10%. This is consistent with other DFT calculations 
that give static polarizabilities in the range 33.1-38.2 for Na2 and 67.1-78.7 for Na4 
(see Ref. and results cited therein). In the case of Na4, the orbital functionals SIC and 
the EXX yield slightly better static polarizabilities than the LDA and the GGAs, but the 
results are still smaller than the experimental value. We also rejoice that the two methods 
used to calculate the static polarizability yield very similar results. The neglect of correlation 
as well as temperature effects is responsible for the obtained smaller polarizabilities values. 
We indicate that the simple argument that a potential with the correct asymptotic will lead 
to more localized charge and therefore lower polarizability does not hold for these small 
systems. These results give support to previous studies on simple metaUic jellium spheres 



B. Silicon 

We have calculated the two simplest hydrogen terminated silicon clusters: silane (SiH4) 
and disilane (Si2H6). These systems pose a much harder challenge than the simple sodium 
clusters (where even jellium calculations within the LDA yield reasonable results), not only 
because of the presence of the p-electrons, but also due to the hydrogen, which is very hard 
to describe by a reasonably soft-core pseudopotential. As expected, the five xc functionals 
we tried yielded quite dissimilar results. 

Our calculations for the IP and HOMO-LUMO gaps for silane and disilane are summa- 
rized in Table UndoubtfuUy, EXX and LB94 gave the best IPs of all five functionals 
tested, almost at the level of the much more involved GW with exciton effects or Monte- 
Carlo calculations. SIC is slightly worse, and LDA and PBE yield, as usual, completely 
irreal IP. Although the IP changes by more than 4 eV going from the LDA to the EXX, the 
HOMO-LUMO gap increases by just around 0.5 eV: The main difference is, once more, a 
nearly rigid shift of the eigenvalues. However, this shift brings the relevant single particle 
transitions below the ionization threshold. The results are also consistent with the finding 



that EXX calculations provide larger gaps for semiconductors ||3^ . 

We will discuss separately the silane and disilane spectra. The silane spectrum (Fig. ^ 
consists of three peaks between 8 and 12 eV (the two peaks derive from a Jahn- Teller 
splitting of the triply degenerate 2t2 — ^ 4s transition |3^), followed by a much broader 
feature at higher energies. The curves obtained with the two traditional functionals, LDA 
and PBE, are quite similar to each other, and the onset for absorption is underestimated 
by around 1 eV. The underestimation of the onset of absorption is a well known deficiency 
of LDA based functionals that is even more dramatic in the case of infinite bulk systems 
where excitonic effects and band-gap renormalization are not properly described by these 
simple functionals [^]. The SIC spectrum is unphysically shifted to lower energies (fact that 
could be anticipated by looking at the low SIC HOMO-LUMO gap), and the first peak is 
split. This second effect is an artifact of SIC, which, in this implementation, spontaneously 
breaks the 3-fold degeneracy of the HOMO state. LB94 and EXX behave quite well: the 
onset for absorption is now correct, and the error in the position of the first three peaks is 
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reduced by a factor of 2 from the LDA/PBE results (see Table 0). van Leeuwen's functional 
performs marginally better in this case than EXX, which overestimates slightly the excitation 
energies. The transitions close to the ionization threshold give rise to a peak in the spectrum 
that in experiments is usually broader than in the calculations, as new decaying channels 
are available that are not included in the present calculations. Although the functionals we 
used do not properly include the damping of the excitations, finite lifetimes can be simulated 
by convoluting the calculated spectrum with some Lorentzian function. 

The disilane LDA, PBE, LB94 and SIC spectra are all very similar (see Fig. ^ and consist 
of five peaks in the 7-12 eV interval, followed by a broader feature at higher energies. All 
these curves suffer from the same deficiency: The separation between the second and third 
peak is too small, so the second peak appears as a shoulder of the third due to the limited 
resolution of the calculations (in Table |V| the second line for Si2H6 describes the position 
of this shoulder). The curve closest to experiment was calculated with the EXX functional 
(note its first peak at exactly the experimental value). SIC yielded a quite unreasonable 
spectrum, consisting of broad peaks at too low energies, and consolidated its position of the 
worst of the five functionals tested (as what concerns the calculation of optical absorption 
spectra, of course). 

For these clusters we also provide the static polarizabilities in Table For silane 
we find, in agreement with the results for sodium clusters, that all the methods lead to 
polarizabilities lower than the experimental values. In spite of the proper asymptotics of the 
EXX functional, LDA yields the best results, which indicates that a subtle error cancellation 
is taking place, and that it is important to have some sort of correlation in the functional. 
This result seems to hold true for these small systems and deserves a more detailed analysis. 

C. Xc potential vs xc kernel 

In order to study the sensitivity of the optical spectra to the xc functional used to 
calculated the ground-state, and to the time-dependent xc potential used to propagate in 
time, we plot, in Fig. ^, the spectrum of silane, for some choices of v^J and v^J (see Eq. (pTsD ). 
In the first plot we take ground states calculated with several xc functionals (EXX, LB94 
and LDA), and propagate these states with the LDA, while in the second we fix the ground 
state (calculated with the EXX functional), and use three different potentials to propagate 
it (LDA, LB94 and EXX). The results clearly indicate that the optical spectrum of SiH4 
is much more sensitive to the choice of the starting state, than to the quality of v^J. This 
numerical evidence is in good agreement with the findings of Petersilka et al. [Q. From the 
second plot of Fig. |^ we see that the ALDA provides the strongest attractive /xc, shifting 
the transitions to lower energies by as much as 0.4 eV, the effect being stronger for the 
high energy transitions. Nevertheless, the role of these functionals needs to be quantified for 
larger systems and, in particular, for extended systems where the form of the initial potential 
is expected to be less important and where correlations are supposed to be the dominant 
contribution to the excitation spectra. This different behavior of finite and infinite systems 
poses a theoretical problem when designing functionals as the two classes of systems are 
mainly sensitive either to the spatial non locality or to the frequency dependence of the xc 
functional. 
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V. CONCLUSION 



We presented a detailed study of the impact of different xc potentials and kernels in 
describing the static and dynamical properties of small clusters. The ionization threshold 
as well as specific single-particle transitions are very sensitive to the specific form of the 
xc potential used to obtain the ground state of the system. The static polarizabilities 
turn out to be rather close to experiment but always smaller than the observed values. 
Surprisingly, the polarizabilities calculated with the LDA or the GGA are better than the 
ones calculated with LB94, SIC or EXX. In the case of the EXX, this is related to the 
absence of dynamical correlations, while the LB94 and the SIC functionals seem to break 
the subtle error cancellation between exchange and correlation always present in the LDA 
and the GGA functionals. 

For the optical spectra a different conclusion is reached. EXX and LB94 perform remark- 
ably well, even if the high lying excitations (below the ionization threshold) appear shifted 
to higher energies. This error is traced back to the weak effective attractive xc kernel. It 
seems clear that the right asymptotics are important but not enough to get a more quan- 
titative agreement with experimental data. SIC performs badly, especially for the silicon 
clusters. The more traditional functionals, the LDA and the GGA, showed their well known 
problems: for example, the onset for optical absorption is underestimated by around 1.5 eV 
in SiH4. 

Finally, our spectra calculated with different xc-functional for the ground state calcu- 
lation and for the time evolution {v^J ^ v^c )^ suggest that the choice xc potential used 
to calculated the ground state is more important than the choice of the xc kernel for the 
case of finite systems. This conclusion should not be valid when one is dealing with infinite 
periodic systems. 
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TABLE CAPTIONS 



Table |: Calculated ionization potentials, (IP), obtained from the highest occupied KS 
eigenvalue, and HOMO-LUMO gaps for Na2 and Na4 for different xc functionals compared 
to available experimental data ||3^. All values are in eV. 

Table |I[ Static polarizabilities for the small sodium clusters estimated using a finite 
electrical field (FF), or calculated from the uj = Fourier transform of the time- dependent 
dipole moment (FT). The results are compared to experiments |^ and other calculations 
113 2||. All values are in A3. 

Table Calculated ionization potentials (IP) obtained from the highest occupied KS 
eigenvalue, and HOMO-LUMO gaps for SiH4 and Si2H6. The available experimental data 
as well as results from many body calculations are given for comparison. DMC stands for 
the diffusion quantum Monte-Carlo calculations of Ref. GW for the GW with exciton 
effects results (GW-Bethe-Salpeter equation) of Ref. and finally, QMC for quantum 
Monte-Carlo calculations EOH. All values are in eV. 



Table |Vf Static dipolar polarizabilities for the hydrogenated silicon clusters obtained 
from the u = Fourier transform of the time-dependent dipole moment. All values are in 

A3. 

Table |V|: Lowest excitation energies for SiH4 and Si2H6 obtained with different xc- 
functionals compared to experiments and other calculations based on a many body for- 
malism. The SIC results are absent due to the difficulty to identify the peaks in the spectra. 
DMC stands for the diffusion quantum Monte-Carlo calculations of Ref. |3^, and the calcu- 
lated value of 9.47 should be compared with the average of the two first, Jahn- Teller splitted, 
levels. GW are the GW with exciton effects results (GW-Bethe-Salpeter equation) of Ref. 



40| , and finally, QMC for Quantum Monte-Carlo calculations [0. All values are in eV 
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FIGURE CAPTIONS 



The "exp" 
The experimental 



Fig. |I|: Averaged dipole strength for Na2, calculated for several xc functionals 
curve is from Ref . , who adapted the experimental results of Ref . [pO[ 
curve is plotted in arbitrary units. 

Fig. 1^: Averaged dipole strength for Na4, calculated for several xc functionals. 

while 



exp 
"GW" is the 



stands for the experimental photodepletion cross sections of Ref. |^ 
many-body calculation of Ref. |^T|, which includes self-energy and excitonic effects. These 
two curves are shown in arbitrary units for the sake of comparison. 

Fig. ^: Averaged dipole strength for SiH4, calculated for several xc functionals. The 
"exp" curve is from Ref. |^ . 

Fig. ^: Averaged dipole strength for Si2H6, calculated for several xc functionals. The 
"exp" curve is from Ref. p4| . 

Fig. ^: Averaged dipole strength for SiH4, calculated using different xc-functionals fro 
the calculation of the ground state and for the time evolution (v^J ^ v^^^ ) (see Eq. |13[). 
For each curve, the first functional mentioned corresponds to v^J and the second to v^^ 
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TABLES 



TABLE L 





PZ 


PBE 


LB94 


SIC 


EXX 


exp 


371 


Na2 IP 

Na2 HOMO-LUMO 


3.21 
1.36 


3.22 
1.36 


4.95 
1.25 


5.17 
1.51 


4.76 
1.45 


4.89 


Na4 IP 

Na4 HOMO-LUMO 


2.75 
0.60 


2.75 
0.63 


4.44 
0.36 


3.85 
0.47 


3.88 
0.72 


4.27 



TABLE II. 





PZ 


PBE 


LB94 


SIC 


EXX 


Ref. [32| 


exp 


Na2 FF 
Na2 FT 


35.0 
34.9 


34.3 
34.1 


31.8 
31.6 


33.7 
33.5 


34.9 
34.7 


37.0 


39.3 
39.3 


Na4 FF 
Na4 FT 


76.7 
75.9 


75.4 
74.6 


73.4 
71.7 


80.0 
78.6 


77.4 
76.3 


78.7 


83.8 
83.8 



TABLE III. 





PZ 


PBE 


LB94 


SIC 


EXX 


DMC [39| 


GW 1 


m 


QMC |4C| 


exp 


SiH4 IP 


8.53 


8.55 


12.7 


11.8 


13.1 


12.88 


12.7 


12.6 


12.61 HI 


SiH4 HOMO-LUMO 


8.10 


8.12 


8.40 


7.70 


8.77 




13.0 






SiaHg IP 


7.40 


7.37 


11.2 


9.95 


10.9 


10.90 








10.53-10.7 m 


Si2H6 HOMO-LUMO 


6.76 


6.80 


6.58 


5.98 


7.17 













TABLE IV. 





PZ 


PBE 


LB94 


SIC 


EXX 


exp 


431 


SiH4 


5.11 


4.94 


4.66 


4.71 


4.53 


5.44 


SisHe 


9.87 


9.58 


9.27 


10.0 


8.70 





TABLE V. 





PZ 


PBE 


LB94 


EXX 


DMC 1391 


GW 1401 


QMC @ 


exp 441 


SiH4 1st 
SiH4 2nd 
SiH4 3rd 


8.23 
9.30 
10.1 


8.25 
9.39 
10.1 


8.75 
9.80 
11.2 


8.93 
9.98 
11.1 


9.47 


9.2 


9.1 


8.8 
9.7 
10.7 


(|A1) 


0.52 


0.49 


0.22 


0.27 










Si2H6 1st 
SizHe 2nd 


7.30 
8.50 


7.35 
8.60 


7.38 
8.42 


7.60 
8.63 








7.6 
8.4 


(|A|) 


0.20 


0.23 


0.12 


0.12 
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FIGURES 



FIG. 1. 
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FIG. 2. 
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FIG. 3. 
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FIG. 4. 
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FIG. 5. 
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